o 



ctf 



CN 

in 

o 
o 



13 



Research in Astron. Astrophys. 2010 Vol. X No. XX, 000-000 D , . 

riesearch m 
http://www.raa-journal.org http://www.iop.org/journals/raa Asfrnnnmv inri 

Astrophysics 



Nonlinear time series anaysis of the light curves from the black 
hole system GRS1915+105 

K. P. Harikrishnan 1 , R. Misra 2 and G. Ambika 3 

1 Department of Physics, The Cochin College, Cochin-682002, India; kp_hk20Q2@yahoo.co.in 

2 Inter University Centre for Astronomy and Astrophysics, Pune-41 1007, India; 



—». , rmtsraWiucaa.ernet.in 

*vi ■ 3 Indian Institute of Science Education and Research, Pune-41 1021, India; g.ambika@iiserputie.ac. in 

Received [year] [month] [day]; accepted [year] [month] [day] 

ta" 

Abstract GRS 1915+105 is a prominent black hole system exhibiting variability over a 
,-Q • wide range of time scales and its observed light curves have been classified into 12 tempo- 

♦pH' ral states. Here we undertake a complete analysis of these light curves from all the states 

Q . using various quantifiers from nonlinear time series analysis, such as, the correlation di- 

mension (-D2), the correlation entropy (K2), singular value decomposition (SVD)and the 
multifractal spectrum (f(a) spectrum). An important aspect of our analysis is that, for es- 
timating these quantifiers, we use algorithmic schemes which we have proposed recently 
and tested successfully on synthetic as well as practical time series from various fields. 
Though the schemes are based on the conventional delay embedding technique, they are 
!> ' automated so that the above quantitative measures can be computed using conditions pre- 

scribed by the algorithm and without any intermediate subjective analysis. We show that 
nearly half of the 12 temporal states exhibit deviation from randomness and their complex 
temporal behavior could be approximated by a few (3 or 4) coupled ordinary nonlinear 
differential equations. These results could be important for a better understanding of the 
processes that generate the light curves and hence for modelling the temporal behavior of 
such complex systems. To our knowledge, this is the first complete analysis of an astro- 
physical object (let alone a black hole system) using various techniques from nonlinear 
dynamics. 
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1 INTRODUCTION 

Most of the systems in Nature are described by models which are inherently nonlinear. Since the discov- 
ery of deterministic chaos a few decades back and the development of various techniques in subsequent 
years, there remained the exciting prospect of a better understanding of the complex behavior shown by 
various natural systems in terms of simple nonlinear models. Evidence for low dimensional chaos has 
been reported - and disputed - not only in physical scienc es, but also in many other fields such as, phys- 
iology, economics and social sciences dSchreiberlll999h . Particular attention has been paid to systems 
producing strange and chaotic attractors, the word strange refering to metric properties such as fractal 
dimension and the word chaotic representing dynamic properties like exponential divergence of nearby 
trajectories in phase space. A large number of techniques and measures from nonlinear dynamics and 
chaos theory are routinely being employed for the analysis of such systems. Excellent text books are 
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now available that give a background knowledge on various methods in nonlinear dynamics dHilbornl 
1 1 994 ISprottl 120031: iLakshmanan & Rajasekaill2003l) . 

One major difficulty in the analysis of real world systems is that our knowledge regarding the 
system is usually limited to a single scalar variable recorded as a function of time, called the time 
series. Therefore, a great deal of effort has been devoted to the characterisation of underlying attractors 
reconstructed from time series. The large number of tech niques and computational schemes used for 
this purpose have be en discussed in detail by many authors JKantz & Schreiben,ll997l : lAberbane]|,ll996t 
lHeggeretallll999h . 

Among the most important quantifiers used for the analysis of time series data are the correlation 
dimension (D2), the correlation entropy (if 2) and the multifractal spectrum. Correlation dimension is 
often used as a discriminating statistic for hypothesis testing to detect nontrivial structures in the time 
series. But when the time series involves colored noise, a better discriminating measure is considered to 
be K2 JKennel & Isabella, 1 19921) . Finally, a complete characterisation of the underlying chaotic attractor 
is done using the generalised dimensions D q and the f(a) spectrum. 

We have recently proposed automated algorithmic schemes dHarikrishnan et all 120061 l2009al) for 
the computation of D2 and K2 from time series based on the delay embedding technique and applied it 
successfully to various types of time series data including that from standard chaotic systems, data 
added with white and colored noise and practcal time series like EEG and ECG. A generalisation 
of these schemes to compute the multifractal spectrum of a chaotic attractor has also been proposed 
dHarikrishnan et alll2010ll2009bl) . These schemes provide a nonsubjective approach for the characteri- 
sation of strange attractors inherent in time series. 

It should be noted that so far, most of the analysis of the light curves from X-ray binaries and active 
galactic nuclei (AGN) have used the conventional techniques such as the power spectrum and distribu- 
tion. It is widely believed that the light intensity variations are mostly stochastic in nature. For example, 
it has been shown in the case of the most prominent black hole system, Cygnus X-l, that the observed 
light curves, at least in certain ti me scales, are cons istent with some static nonlinear transformatons of 
stochastic variations in intensity dUttlev et all 120051) . The authors argue that models based on nonlinear 
dynamics are not required to explain the data. 

But there are also some analysis based on nonlinearity measu res that have been attempted earlier 
dVoses et all Il987t iNorris & Matilskvl Il989t iTimmeret all |2000|) on the light curves of some promi- 
nent black hole systems, such as, Her X-l and Cygnus X-l. But these studies have so far not been able 
to provide conclusive evidence for nontrivial structures in the temporal behavior of such systems. One 
reason for this has been the limited number o f data sets available from such sources with sufficient 
signal to noise ratio required for such analysis dNorris & Matilskvl 119891) . The scenario has changed in 
the last few years as enough data are now available through RXTE observations. Recently, a nonlinear 
time series analysis performed on light intensity data from the white dwarf variable PG 135 1+489 has 
enabled much more i nform ation regarding the system compared to the conventional power spectrum 
analysis dJevtic et alll2005l) a nd attempts have be en made to use these analysis to differentiate between 
neutron stars and black holes dKarak et all 1201 Ol) 

Studies on GRS1915+105 have been limited because it became active only over a decade ago. But 
the system turns out to be unique among all such sources in that it seems to flip from one state to another 
continuously with each state having its own temporal variability different from the other st ates. The light 
curves have been classified into 12 spectroscopic classes based on RXTE observations bv lBelloni et al.l 
d2000l) . The nature of the light curves changes completely as the system flips from one state to another. 
Evidently, pure stochastic processes cannot account for such qualitative changes in the light curves. 
Hence the question naturally arises whether some nonlinear deterministic processes are also involved. 
We do find evidence to support this. 

In this paper, we apply the above mentioned automated schemes developed by us to undertake a 
complete analysis of the X-ray light curves from GRS 1 9 1 5+ 1 05 . Earlier, a surrogate analysis with D^ as 
the discriminating measure has shown that a few of these states manifest the t ime evolutions ana logous 
to that from low dimensional nonlinear systems with some inherent noise dMisra et all 120061) . This 
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Fig. 1 The D^ values of the Rossler attractor (with error bar), as a function of M along with 
the Z?2 values of 20 surrogates (dashed lines). All computations are done with 10000 data 
points. 



motivates us to undertake an exhaustive numerical analysis of the light curves from the source in all its 
temporal states using the prominent tools from nonlinear dynamics. 

Another motivation for the present investigation has been derived from the fact that the accretion 
disk in such systems are driven by magneto-hydrodynamic turbulence which is an intrinsically non- 
linear process. A model for such a process should be nonlinear and is expected to show qualitative 
changes in its behavior as a control parameter is varied. For the X-ray radiation from an accretion disk, 
the rate of mass accretion could poss ibly be considered as a suitable cont rol parameter. There also ex- 
ists theoretical models to this effect dVoges et all 1 19871 : lAtmanspacher et all Il989al) from which it is 
possible to derive the temporal variability of the X-ray radiation in different regimes of mass accretion 
rate. Since many of the states of the black hole system under study show nonlinear character, a com- 
plete analysis of the light curves using various nonlinear measures can greatly help in the search for a 
nonlinear deterministic model to describe the temporal variability of the system. 

Our paper is organised as follows: All the quantitative measures used in this paper and the corre- 
sponding computational schemes are discussed in detail in the following section. While §2.1 and §2.2 
present the computational details for D2 and K2, §2.3 and 2.4 are concerned with SVD and f(a) spec- 
trum respectively. The time series from a standard chaotic system - the Rossler system - is used as an 
example to illustrate the results in all the cases. The analysis of the X-ray light curves from the black 
hole system is then undertaken in §3 and the conclusions are drawn in §4. 
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2 QUANTITATIVE MEASURES USED FOR THE ANALYSIS 

2.1 Correlation dimension and surrogate analysis 

Correlation dimension D 2 is often used as a discriminating statistic for hypothesis testing. T he con- 
ventio nal method for the calcu lation of D 2 is the delay embed ding method first introduced by iTakensI 
dl98l[) . and used effecti vely by Grassberger & Procaccial dl983l) . now known as the GP algorithm. More 



details can be found in ISauer et al.l d 199 lb . It creats an embedding space of dimension M with delay 
vectors constructed by splitting a discretely sampled scalar time series s(ti) with delay time r as 

Xi = [s(U), s(U + r), ...., s{U + {M - l) r )] (1) 

The correlation sum is the relative number of points within a distance R from a particular (i th ) data 
point, 

1 N v 

Pi (R)= lim — J2 H(R-\ Xi - Xj \) (2) 

where N v is the total number of reconstructed vectors and H is the Heaviside step function. Averaging 
this quantity over N c randomly selected x^ or centers gives the correlation function 

1 Nc 
C M {R) = —Y.P^ R ) P) 

i 

The correlation dimension D 2 (M) is then defined to be, 

D 2 = lim d(logC M (R))/d(log{R)) (4) 

which is the scaling index of the variation of Cm(R) with R as R — > 0. In pratice, a linear part in the 
logCniR) versus logR plot is identified subjectively, called the scaling region, and its slope is taken 
as D 2 . But in our computa tional scheme, this is done algorithmically as discussed in detail elsewhere 
( iHarikrishnan et aLld2006l) ) and the scheme computes D 2 with error bar as a function of M. The scheme 
has also been shown to be suitable for hypothesis testing using surrogate data. 

The rationale behind surrogate analysis is to formulate a null hypothesis that the data has been 
generated by a stationary linear stochastic process, and then attempt to reject it by comparing a suit- 
able measure for the data with appropriate realisations of surrogate data. T he method fo r the gener- 
ation of surrogate data was originally proposed by Theiler and co workers dTheiler et all 1 19921) with 
the Amplitude Adjusted Fourier Transform (A AFT) algorithm. But lSchreiber & Schmitd d 19961 l2000h 
have proposed another iterative scheme, known as IAAFT scheme, which is similar but reported to be 
more consistent in representing null hypothesis dKugiumtzisl 1 19991) for a wide class of stochastic pro- 



cesses. In this work , we apply this scheme to generate surrogate data sets using the TISEAN package 
dHeggeretal.1119991) . 

We first apply the D 2 analysis on different types of data sets such as, standard chaotic time series, 
pure noise and chaotic data added with noise. The time series from the standard Rossler attractor, with 
parameter values a = b = 0.2 and c = 7.8, is used as a reference to test all the computational schemes 
presented in this work. All computations are done with 10000 long data points and 20 surrogates for 
each data. In Fig.l, D 2 of Rossler attractor and surrogates are computed as a function of the embedding 
dimension M, where as in Fig. 2, the same is shown for two pure colored noise data sets with spectral 
index s = 1.5 and 2.0. As expected, the Rossler data show clear deviation from the surrogates while for 
the latter, the null hypothesis cannot be rejected. 

Now the real world data is often contaminated with noise and the question that arises naturally is 
how much amount of noise can suppress the nonlinear component that may be present in the time series. 
In order to study the effect of noise on D 2 using our scheme, we generate two data sets by adding 50% 
of white and coloured noise (with ,s = 2.0) to the time series from the Rossler attractor. The result of 
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Fig.2 The upper panel shows the Z?2 values of pure colored noise with spectral index s = 1.5 
as a function of M along with surrogates. The lower panel shows the same for s = 2.0. 



applying our scheme on these data and their surrogates is shown in Fig. 3. It is found that when white 
noise is added to the system, D 2 of the data increases and for a contamination level of > 20%, it is 
difficult to distinguish between the data and the surrogates. But for colored noise contamination, the 
data is distinguishable from the surrogates for an added noise level of upto 50%. Note that 50% noise 
here means that the noise amplitude is half of that of data. The above results indicate that from a D 2 
analysis, it is difficult to distinguish even moderate amount of noise contamination in a chaotic data. A 
better quantitative measure in such a situation is K 2 discussed in the next section. 

In order to get a quantification of the differences in the discriminating measure between the data 
and the surrogates, we u se the normalised mean sigma deviation (nmsd), proposed by us recently 
JHarikrishnan et aU 12006). For D 2 , this is computed using the expression 



nmsd 2 



M„ 



M„ 



M=2 



D 2 {M)- 



< Df""' 



(M)> 



T sur 
'SD 



r {M) 



(5) 



where M max is the maximum embedding dimension for which the analysis is undertaken, < 
D s 2 urr (M) > is the average of D s 2 urr {M) and a s s u £ r (M) is the standard deviation of D s 2 urr (M). We 
have earlier shown that a value of nmsd < 3.0 implies eithe r white or colored noise d omination in the 
data and the null hypothesis cannot be rejected (see for e.g. lHarikrishnan et all 120061) . It is found that 
for the Rossler attractor data shown in Fig.l, the nmsd = 36.1 and for the two pure colored noise in 
Fig.2, the nmsd = 0.68 and 2.0 for s = 1.5 and 2.0 respectively. For data contaminated by noise in Fig. 3, 
the values are 1.8 for white noise and 3.5 for colored noise. 
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Fig. 3 The figure shows the D 2 values as a function of M for the Rossler attractor data added 
with 50% white noise (s = 0.0) and the same percentage of colored noise (s = 2.0) along 
with their respective surrogates. Note that, in the case of colored noise contamination (lower 
panel), the data values are still below that of the surrogates. 



2.2 Correlation entropy 

The use of K 2 has been limited compared to D 2 for the analysis of time series data. But in cases 
whe re time series involv e colored noise, K 2 is a more effective discriminating measure compared to 
D 2 JRedaelli et all 12002). While D 2 is a geometric measure of the underlying chaotic attractor, K 2 is a 
dynamic measu re represen ting the rate at which information needs to be created as the chaotic attractor 
evolves in time (lOttLll993l) . The standard method for the computation of K 2 is also the delay embedding 
technique. Since K 2 measures the rate at which the trajectory segments are increased as M increases, it 
can be related to the correlation sum Cm (R) by the expression 



Cm{R) oc e 



-MK 2 At 



(6) 



where At is the time step between successive values in the time series. From above, a formal expression 
for K 2 can be written as, 



K 2 At=lim lim lim (-logC M (R)/M) 

H->0 M ->oo JV->oo 



Alternately, K 2 can also be obtained as 
K 2 At 



lim lim lim log(C M {R)/C M +i{R)) 

R— >0 M — >oo AT— >oo 



(7) 



(8) 



Our no nsubjective scheme has been extended for the computation of K 2 as well JHarikrishnan et aU 
l2009al) and we apply that scheme in this work. 
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Fig. 4 The X2 values of the Rossler attractor as a function of M computed using our scheme. 
The values are computed per sec and converges very close to the standard value 1.04 ± 0.08. 



Fig. 4 shows Ki for the Rossler attractor as a function of M computed from the time series using 
our scheme. To show the effect of noise on K%, we generate four different time series by adding 50% 
and 100% white as well as colored noise to the Rossler attractor data. The result of applying our scheme 
to these data sets is shown in Fig. 5. It is clear that while the saturated K-% value increases with the white 
noise addition, the effect of colored noise is quite the opposite. With the increase in colored noise, the 
saturated value of K2 —¥ 0. 

Our scheme can be used for surrogate analysis with K 2 as the discriminating measure as well 
and nmsd can be comput ed using a similar expression as Eq. [5] We have applied this recently 
dHarikrishnan et all l2009al) to the Rossler attractor data with different percentages of white and col- 
ored noise added. For 50% white noise contamination, nmsd with K% as the discriminating measure is 
found to be 4.3, while for the same percentage of colored noise, the value is 2.2. Thus, while the white 
noise contamination can be easily identified through D2 analysis, the presence of colored noise can be 
better inferred by computing K%. 



2.3 Singular value decomposition 

The singular value decom position (SVD) is another i mportant technique used in nonlinear time series 
analys is, first proposed bv lBroomhead & Kind ( 119 86') and for a recent review, see lAthanasiu & Pavlosl 
J200ll) . The method makes use of a trajectory matrix constructed from the experimental time series with 
the rows of the matrix constituting the state vectors in the embedding space. It is then diagonalised 
to find the dominant eigen values and eigen vectors which are used to represented the dynamics. The 
number of dominant eigen values determine the minimum number of dimensions required to unfold 
the complete dynamics and the corresponding eigen vectors give the projections of the reconstructed 
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Fig. 5 Variation of K 2 with M for data sets obtained by adding different amounts of white 
and colored noise to Rossler attractor data. The upper panel shows the result of addition of 
(a) 50% and (6)100% white noise, while the lower panel shows the results with (a)50% and 
(b) 100% colored noise with s = 2.0. It is clear that while K% increases with addition of white 
noise, Ki — > as the percentage of colored noise increases. 



attractors. With such a SVD projection (also called BK projection), one can vis ualise the qua l itative 
nature of the reconstructed attractors. Here we use the standard TISEAN algorithm dHegger et all 1 19991) 
for the computation of BK projections. 

For example, the SVD projection for the Rossler attractor is shown in Fig. 6 (upper panel). To show 
the effect of colored noise, on the SVD projection, we also show in Fig. 6 (lower panel) the BK projection 
for the Rossler attractor added with 50% colored noise with spectral index s = 2.0. It is evident that 
even such large amount of colored noise does not destroy the attractor completely. 



2.4 Multifractal spectrum 

The interest of the multifractal formalism in connection with dynamical systems rests on the fact that 
it provides us with a very efficient method to determine the existence of strange attractors and allows 
a statistical description of these sets. A strange or chaotic attractor normally possess a multifractal 
structure as a result of the stretching and folding of the trajectories in dif ferent directions in phas e space . 
Hence they can be characterised by a spectrum of dimensions D(q) dHentschel & Procaccial Il983l) . 
where the index q can vary from —00 to +00. The clustered regions of the attractor are characterised by 
D(q) values with q > and rarefied regions by D(q) with q < 0, with D(0) giving the simple fractal 
dimension of the set. 

A more convenient method to represent the global scaling properties of the attractor is by using a 
spectrum of singularities characterising the probability measure on the strange attractor. For this, one 
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Fig. 6 The SVD plot of the pure Rossler attractor (upper panel) along with a plot of the 
attractor obtained by adding 50% colored noise (s = 2.0) to the Rossler attractor time series. 



considers the number N a (e) of boxes (with edge length e) required to cover the attractor within a small 
range of a given by a <G [a, a + da]. The parameter a is a continuous variable characterising the local 
scaling properties of the fractal set. Its meaning is that of a local scaling exponent: 



Qj = lim logpi(e)/loge 



(9) 



with Pi(e) representing the probability measure. Thus a measures how fast the number of points within 
a box decreases as e is reduced. It therefore measures the strength of a singularity for e — > 0. Assuming 
that a scaling hypothesis holds, one can write 



N a (e)(xe- fia '> 



(10) 



where the exponent f(a) represents the fractal dimension of subsets of strength a. The graph of f(a) 
as a function of a is called the f(a) spectrum which characterises the global scaling properties of the 
fractal set as a function of the local scaling exponents a. The transformation fro m D(q) to f(a) can 
be show n to be a Legendre transformation, for details, see lHalsey et al.l dl986l) and lAtmanspacher et al.l 
(Il989bh . 

To compute the f(a) spectrum from a time series, we first consider the generalised correlation sum 
given by 

1 Nc 1 N v 

c m( r ) = ^wY,(w E H{R-\xi-*i\)) q - l ] l/{q - l) (ID 



i=i,i#j 
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Fig. 7 The upper panel shows the D q spectrum of the Rossler attractor with error bar, com- 
puted from the time series of 10000 data points. To show the accuracy of fitting, the D q 
values (points) are again shown in the lower panel without error bar along with the best fit 
curve (continuous line) computed using our numerical scheme. 



where the Heaviside function H counts how many pairs of points at (xi,Xj) are situated within a 
distance R. The spectrum of dimensions are then determined by the relation 



D(q) = Vimd(logC q M (R))/d(log(R)) 



R^0 



(12) 



The average value of D q with error bar is then calculated from the scaling region by taking different 
values of R, by extending the numerical procedure discussed above for computing Z?2- 

In order to determine the f(a) spectrum, we make use of the computat i onal sc heme recently pro- 
posed by us and applied to several practical time series JHarikrishnan et all l2009bl) . The scheme uses 
an analytical fit for the f(a) function (involving a set of independant parameters) and calculates the 
corresponding D q curve using the Legendre transformation equations. This curve is then fitted to the 
spectrum of D q values computed from the time series. The best fit curve is found by changing the pa- 
rameters of the /(a) fit, which is then u sed to compute the f i nal f(a ) spectrum. The algorithmic details 
of the scheme are presented elsewhere JHarikrishnan et afl, l2009bl) . The multifractal approach has re- 
cently been employed in the analy sis of several practic al time series, an example being the temporal 
variations of the geomagnetic field dHongre et al!ll999l) . 

To illustrate our scheme, it is first used to compute the f(a) spectrum of the Rossler attractor. 
The spectrum of generalised dimensions D q is computed from the time series taking the embedding 
dimension M = 3. Attempting to compute the f(a) spectrum directly from the D q values leads to 
an incomplete f(a) spectrum. This is mainly due to the fact that the errors in the calculation of D q 
makes the Legendre transformation numerically impractical because of the reversal of slopes. Hence 
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Fig. 8 The / (a) spectrum of the Rossler attractor computed from the best fit D q curve shown 
in the previous figure. 



our scheme uses a different procedure. The f(a) function is a single valued function defined between 
the limits of a m in and a m ax- Since the derivative f'(a) = df(a)/da = q is also single valued, it 
follows that f(a) has a single extremum (i.e. a maximum). Moreover, f(a m in) = f(ohnax) = and 
f'{&min) an d f'(ctmax) tend to oo and — oo respectively. A simple function which can satisfy all the 
above necessary conditions is, 



f(a)=A(a 



,) 71 K 



a) 



72 



(13) 



where A, 71, 72, a m in and a max are a set of parameters characterising a particular f(a) curve. The 
D q curve can be computed from this f(a) fit using the inverse Legendre transformation equations for 
a given set of parameters. It is then fitted to the D q spectrum computed from the time series. The 
statistically best fit D q curve is found by adjusting the parameters of the f(a) function, which is then 
used to compute the final f(a) spectrum. The D q spectrum and its best fit curve for the Rossler attractor 
are shown in Fig. 7 and the f(a) spectrum computed from the best fit curve is shown in Fig. 8. Having 
discussed the various measures and schemes for computing them, we now turn to the analysis of the 
black hole system GRS 19 15+105. 



3 ANALYSIS OF THE BLACK HOLE SYSTEM GRS1915+105 

In this section, we apply all the techniques discussed above to analyse the X-ray light curves from the 
black hole binary GRS 19 15+ 105 . The temporal prop erties of the system have been classified into 12 
different spectroscopic classes bv lBelloni et al.l J2000I) based on the RXTE data. Here we have chosen 
representative data sets for each class. The light-curve for an observation was obtained from the standard 
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Fig. 9 Light curves from the twelve temporal states of the black hole system GRS 1 9 1 5+ 1 05 . 
Only a part of the generated light curve is shown for clarity. 



product^] which provide a 0. 125 sees time resolution summed over all energy channels. Standard prod- 
uct light curves have been generated using a pipeline which considers standard filtering criteria and use 
data from the instruments that were reliably on. While standard products may not have optimal spectral 
response matrix or background model, they are more than adequate for light curve analysis, especially 
for bright sources, like GRS 1915+105, where the background is not important. 

The analysis requires continuous data without gaps. For each class, we have extracted two sets of 
continuous segments for the analysis. The light curves have been generated after rebinning to a time 
resolution of 0.5 sees resulting in ~ 5200 to 6400 continuous data points for each segment. Light curves 
with finer time resolutions are more Poisson noise dominated, while larger binning gives less data points. 
Table 1 gives the observation ID, class, number of data points etc. of all the light curves used for the 
analysis. In the last column, we also provide the temporal behavior of the light curve as resulted from 
our analysis. More detail s regarding t he dat a, such as, average count, expected Poisson noise variation, 
etc. are given elsewhere dMisra et all 120041) . 

Fig. 9 shows all the 12 light curves used for the analysis, which are labelled by 12 different symbols 
representing the 12 temporal states of the black hole system. We show only one set of light curve since 
the second set looks identical for all th e states. The s ystem appears to flip from one state to another 
randomly in time. The classification of iBelloni et alj d2000l) is based on a detailed analysis of all the 
light curves from RXTE data using various linear tools. But it is difficult to differentiate the subtle 
temporal features between the light curves with the help of the linear tools, such as, the power spectrum. 



http://heasarc. gsfc.nasa.gov/docs/xte/recipes/stdprod_ guide, html 
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Fig. 10 The power spectrum for the X-ray light curves from GRS 1915+105 in four represen- 
tative states. 



For example, in Fig. 10, we show the power spectrum for four representative states, whose temporal 
properties are different as per our analysis (see Table 1). While /3 and v are candidates for deterministic 
nonlinear behavior, k is possibly a mixture of nonlinearity and colored noise and 7 is purely stochastic. 
But these distinctions are barely evident from the power spectral variations, though the 7 state appears 
more like a white noise, in agreement with our results. This, once again, emphasizes the importance of 
methods based on nonlinear time series analysis for a better understanding of the temporal properties of 
the light curves. 

Recently, we applied surrogate analysis on all these lig ht curves and show ed that more than half 
of these 12 states deviated from a purely stochastic behavior dMisra et al.l 120061) . Here we combine the 
results of computations of D2, K2 and S VD analysis to get a better understanding regarding the nature of 
these light curves. We have done the surrogate analysis with D2 and K2 and the SVD analysis separately 
for the two sets of light curves. But here we only show plots from representative light curves from one 
set as the plots from the second set are similar and the results are qualitatively the same. Figs. 1 1 and 12 
show the results of surrogate analysis, with D 2 as discriminating measure, on eight different states. Of 
the states shown in these figures, it clear that only two states - 7 and <fr - show purely stochastic behavior. 

Fig. 13 presents the results of computation of K2 for six of the above eight states. Since colored 
noise is also expected in the black hole data, surrogate analysis has been performed with K2 as the 
discriminating measure on all the 12 states. The results are shown in Fig. 14 for four of these states. 
While the behavior of /3, and 7 are consistent with earlier analysis, the behavior of k suggests that it 
is contaminated by colored noise. 

To get a quantitative measure, we now compute nmsd using both D2 and A'2 as discriminating 
measures for all the 12 states from the two sets of light curves and the results are shown in Table 2. A 
careful inspection of the Table reveals the following results. The values of nmsd in the two cases suggest 
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Fig. 11 Surrogate analysis with D 2 as a discriminating measure for the light curves from four 
states of GRS 19 15+ 105. Note that only the 7 state is consistent with noise. 



that the temporal properties of the light curves in the two sets are almost identical. Out of the 12 states, 
four (6, 7, 0, x) are completely stochastic or white noise. Of the remaining eight states, three - k, A and 
p - are contaminated by colored noise and the rest - (3, 8, a, v and p - show signatures of deterministic 
nonlinear behavior in their temporal variations. It is generally expected that all the states contain some 
amount of white noise whose percentage may vary. For example, in the case of the state a, the saturated 
values of D 2 and K 2 improves significantly as the resolution time is increased from 0.5 sec to 1 sec. 
This clearly indicates the presence of Poisson white noise in the data. But the surrogate analysis with 
both D 2 and K 2 confirms that the null hypothesis can be rejected for the light curve in the a state. 

We next perform a SVD analysis on all the states which clearly show the qualitative nature of the 
underlying attractors. The plot of attractors for selected states are shown in Fig. 15. The most interesting 
plot is for the p state which shows a typical limit cycle type attractor. Also, note that the SVD plot for 
k has a nontrivial appearance, even though surrogate analysis suggested the presence of colored noise. 
This is also true in the case of the other two identical states, A and p. Thus these three states could be a 
mixture of deterministic nonlinearity and colored noise. Thus, based on our results, the 12 states can be 
divided into 3 broader classes from the point of view of their temporal properties. It turns out that some 
of these states which are spectroscopically different, behave identically in their nonlinear dynamics 
characteristics. This may be an indication of some common features in the mechanism of production 
of light curves from these states. The temporal behavior of each state as obtained from our analysis is 
indicated in the last column in Table 1 . Since the behavior is identical for two different observation IDs 
in all cases, it may be concluded that the results pres ented here are not de pendent on sample selection 
and are applicable for all the light curves classified bv lBelloni et al.l d2000h . 

Finally, we show the results of multifractal analysis of all the light curves except the four which 
show purely stochastic behavior and hence the f(a) spectrum is irrelevant. Our non subjective scheme 
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Fig. 12 Same as the previous figure, but with four other states. Again, only one state 
shows random behavoir. 



for computing D q and /(a) spectrum discussed in Sec. 2.4, provides us with a set of parameters that 
can be used to compare the fractal properties between different states as reflected in the light curves. To 
show the details of the computations, we first take a typical state. In Fig. 16, we show the D q spectrum 
for the v state along with the best fit curve from which the /(a) is computed as given in Fig. 17. This 
is repeated for the other states as well and the results for four other states are shown in Fig. 18. The 
multifractal nature of the attractors is evident from the figures. Since the computation is done under 
fixed conditions prescribed by the algorithmic scheme, the associated parameters characterising the 
spectra can give a better representation for comparison between various states. Our results indicate that 
the spectra and the associated parameters are typically different for each state and do not show any clear 
trend among members that display strong deterministic nonlinear behavior. This can be inferred from 
Fig. 18 for the case of range of scales. Thus, it turns out that there are subtle differences between the 
states belonging to the same dynamic class with respect to multifractal scaling as well, apart from linear 
spectral characteristics based on which the 12 states are divided. 



4 DISCUSSION AND CONCLUSION 

Identifying nontrivial structures in real world systems is considered to be a challenging task as it re- 
quires a succession of tests using various quantitative measures. Eventhough a large number of potential 
systems from various fields have been analysed so far, the results remain inconclusive in many cases. 
Here we present an example of a very interesting astrophysical system, which we analyse using several 
important quantifiers of low dimensional chaos. By using the time series from a standard chaotic system 
- the Rossler attractor - we first test the computational schemes used for the analysis. These schemes 
are then applied to the light curves from the black hole system. We find that out of the 12 spectroscopic 
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Fig. 13 Variation of K 2 as a function of M for the light curves from six different states of 
GRS1915+105. While the K 2 values of four states converge much like a low dimensional 
chaotic system, K 2 for the k state continue to decrease as M increases indicating colored 
noise contamination. Though K 2 for 7 state converge, its value is much higher and close to 
that of the white noise. 



states, only four are purely stochastic. The remaining show signatures of deterministic nonlinearity, 
with three of them contaminated by colored noise. All these 8 states are found to have D 2 < 4 so that 
their complex temporal behavior can be approximated by 3 or 4 coupled ordinary differential equations. 
Based on our results, the 12 states can be broadly classified into three from a dynamical perspective: 
purely stochastic with D 2 — > 00, affected by colored noise and those which are potential candidates for 
low dimensional (D 2 < 4) chaotic behavior. 

It should be noted that lBelloni et al.l J2000I) classified the light curves into 12 states based on their 
count rate, variability and spectral characteristics. In other words, this is a classifi cation based on linea r 
characteristics of the light curves. Ours is not a classification as in the strict sense of lBelloni et al.ld2000l) . 
We only show that some of the light curves which appear different based on their variability and spectral 
properties can be grouped together when viewed from a dynamical perspective. 

Our results could be significant in many ways. First of all, this is the first real evidence of a pos- 
sible multifractal attractor in the time series of a black hole system. The fact that the light curves from 
many of the temporal states have underlying strange attractor like behavior, increases the possibility 
that the temporal variability in the time scales within these states are governed by some inherent nonlin- 
ear processes with a few degrees of freedom. In other words, the complex nonlinear partial differential 
equations that are known to govern the hydrodynamic flow can be approximated by a set of ordinary 
differential equations and hence can be more easily studied and understood. Moreover, the result that 
some of the states which are spectroscopically different, have approximately the same nonlinear char- 
acteristics is interesting from certain dynamical aspects, such as the mechanism of production of light 
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Fig. 14 Surrogate analysis of the light curves from four states of the black hole system with 
K2 as the discriminating measure. Note that while data and the surrogates can be distin- 
guished for /3 and 6, k and 7 behaves much like colored noise and white noise respectively. 



curves. It is well known that GRS1915+105 is a unique black hole system with many temporal states 
which vary over a wide range of time scales. Many of the questions regarding this variability and the 
exact mechanism of production of light curves still remain unanswered. 

Another question is regarding the structure of variability between the 12 spectroscopic states. It 
has been suggested that all the observed light curves could be interpreted in term s of thr ee basic states 
(a hard state and two softer states) and a sequence of transitions between them (Belloniet a l.. 2000). 
This could, in principle, give rise to a much larger variety of light curves. But the system chooses only a 
handful of these sequences. This possibly suggests that the structure of time variability is not random, but 
controlled by some physical parameter which must be connected to the basic properties of the accretion 
disk. The presence of deterministic nonlinear behavior in the system further substantiates this idea. 

Interestingly, the system appears to be in the \ state for most of the time which is identified as 
purely stochastic in our analysis. We have also analysed many samples of the long time average of the 
light curves from the system and found that they show purely random bahavior. Thus, one possibility 
is that the states other than the \ state ma Y we H be short time flips due to some changes taking place 
within the system. But many of these short time states pick up much less amount of noise revealing, for 
example, the underlying nonlinear character. Thus, an interesting question is whether the states such as 
/3 and 6 have a different underlying mechanism of production of light curves compared to the x state. 
Or, is the excessive amount of white noise in that state itself is suppressing the nonlinear properties? 
The question whether the different states are temporal manifestations of a single underlying mechanism 
or are they dynamically different, will be vital for a proper modelling of this fascinating astrophysical 
object. 
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Fig. 15 The plot of the attractors underlying four states of the black hole system reconstructed 
via SVD analysis. Except the <j> state, which behaves similar to a white noise, all the others 
indicate the presence of underlying attractors, the most interesting being the p state. 
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Table 1 Details of the light curves from GRS 1915+105, in all the 12 
spectral classes, used for the analysis. For each class, light curves from two 
Observation IDs have been analysed, as indicated. The second column gives 
the number of continuous data points after rebinning. The last column indi- 
cates the temporal behavior of the class as obtained from our analysis. 



Obs. ID 



No. of Data Points Class 



Temporal Bahavior 



10408-01-10-00 
20402-01-46-00 


6146 
5504 


20402-01-45-00 
10408-01-15-00 


6156 
5764 


20402-01-03-00 
20402-01-31-00 


6244 
5876 


10408-01-40-00 
10408-01-41-00 


6024 
5312 


20187-02-01-00 
20402-01-22-00 


6010 

5220 


20402-01-33-00 
20402-01-35-00 


6240 
6244 


20402-01-37-00 
20402-01-36-00 


6080 
5648 


10408-01-08-00 
10408-01-34-00 


5688 
5756 


10408-01-17-00 
20402-01-41-00 


6010 
5466 


20402-01-56-00 
20402-01-39-00 


6324 
6180 


10408-01-12-00 
10408-01-09-00 


6286 
5580 


10408-01-22-00 
20402-01-04-00 


6022 
5382 



Deterministic Nonlinear 



Deterministic Nonlinear 



Deterministic Nonlinear 



Deterministic Nonlinear 



Deterministic Nonlinear 



k Deterministic Nonlinear + Colored Noise 



A Deterministic Nonlinear + Colored Noise 



/i Deterministic Nonlinear + Colored Noise 



White Noise 



White Noise 



White Noise 



White Noise 



Timmer, J., Schwarz, U., Voss, H. U., Wardinski, I., Belloni, T., Hasinger, G., van der Klis, M. & Kurths, J. 2000, 

Phys. Rev. E, 61, 1342. 
Jevtic, N., Zelechoski, S., Feldman, H., Peterson, C. & Schweitzer, J. S. 2005, Astrophys. J., 635, 527. 
Karak, B. B., Dutta, J., & Mukhopadhyay, B. 2010, Astrophys. J., 708, 862 
Belloni, T, Klein- Wolt, M., Mendez, M., van der Klis, M. & van Paradjis, J. 2000, Astron. and Astrophys., 355, 

271. 
Misra, R., Harikrishnan, K. P., Ambika, G.& Kembhavi, A. K. 2006, Astrophys. J., 643, 1 1 14. 
Atmanspacher, H., Demmel, V., Morfill, G., Scheingraber, H., Voges W. & Wiedenmann, G., in Measures of 

Complexity and Chaos, Eds. Abraham, N. B., Albano, A. M., Passamante, A. & Rapp, R. E. 1989a, (Plenum 

Press, New York). 
Takens, F. 1981, Lecture Notes in Mathematics, Vol. 898, (Springer, New York). 
Grassberger, P. & Procaccia, I. 1983, Phys. Rev. Lett., 50, 346; Physica D, 9, 189. 
Sauer, T, Yorke, J. & Casdagli, M. 1991, /. Stat. Phys., 65, 579. 

Theiler, J., Eubank, S., Longtin, A., Galdrikian, B. & Farmer, J. D. 1992, Physica D, 58, 77. 
Schreiber, T. & Schmitz, A. 1996, Phys. Rev. Lett., 77, 635. 



20 



K. P. Harikrishnan, R. Misra & G. Ambika 



a 




Fig. 16 The spectrum of generalised dimensions with error bar (upper panel) for the v state 
corresponding to embedding dimension M = 3. The lower panel shows the D q values without 
error bar and the best fit curve. 



Table 2 The values of nmsd computed for all the 12 GRS states from the 
two sets of light curves with D-2 and K% as discriminating measures. 
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Fig. 17 The multifractal spectrum for the v state computed from the best fit curve for D q 
shown in the previous figure. 
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Fig. 18 The /(a) spectrum for the light curves corresponding to four states of the black hole 
system computed using our scheme with M = 3. Note that GRS state has been labelled alpha 
in order to avoid confusion with the scaling index a. 



